libname mort "E:\research\c.NE_mort_temp";

* the final model;
proc phreg data=mort; 
class sex01 race3 agec5 time;
model count*dead(0)= sumtemp wintemp sdsumtemp sdwintemp pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

*without mutual adjustment;
proc phreg data=mort; 
class sex01 race3 agec5 time;
model count*dead(0)= sumtemp wintemp pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

proc phreg data=mort; 
class sex01 race3 agec5 time;
model count*dead(0)= sdsumtemp sdwintemp pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

***********************************************************************
***The percentiles distribution of the purely spatial variation**********;
proc summary data=mort.mort_temp_ses_smk_all;
class zipc;
var sumtemp wintemp sdsumtemp sdwintemp;
output out=new mean= ;
run;

data mort.temp_spatial;
set new;
if zipc ne .;
drop _type_ _freq_;
run;

proc univariate data=mort.temp_spatial;
var sumtemp wintemp sdsumtemp sdwintemp;
run;

***The percentiles distribution of the purely temporal variation**********;
proc glm data=mort.mort_temp_ses_smk_all;
class zipc;
model sumtemp=zipc;
output out=ressum r=rsum;
run;

data spatial;
set mort.temp_spatial(rename=(sumtemp=zipmean_sum wintemp=zipmean_win sdsumtemp=zipmean_sdsum sdwintemp=zipmean_sdwin));
run;

proc sort data=spatial;
by zipc;
run;

proc sort data=mort_fips_smk;
by zipc;
run;

data mort2;
merge mort_fips_smk (in=a) spatial(in=b);
by zipc;
if a;
run;

data mort_final;
set mort2;
sumdif=sumtemp-zipmean_sum;
windif=wintemp-zipmean_win;
sdsumdif=sdsumtemp-zipmean_sdsum;
sdwindif=sdwintemp-zipmean_sdwin;
run;

data mort.mort_temp_ses_smk_dif;
set mort_final;
run;

proc phreg data=mort_final; 
class sex01 race3 agec5 time;
model count*dead(0)= zipmean_sum zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

*distribution;
data distri;
set mort;
keep zipc year sumtemp wintemp sdsumtemp sdwintemp zipmean_sum zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif;
run;
proc sort data=distri nodupkey;
by zipc year;
run;
proc means data=distri mean median std min max range qrange q1 q3;
run;

*correlation matrix;
data corr;
set distri;
if sumdif ne .;
run;
proc corr data=corr;
var sumtemp wintemp sdsumtemp sdwintemp zipmean_sum zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif;
run;


*testing interactions;
*Repeat these models for each exposure index;
data new.mort_temp_ses_all; 
set new.mort_temp_ses_all; 
age01=1;
if age<75 then age01=0;
run;

proc phreg data=new.mort_temp_ses_all; 
class sex01 race3 agec5 time age01(ref='1') /param=ref;*repeat for each reference group;
model count*dead(0)= sumtemp wintemp sdsumtemp sdwintemp sumtemp*age01 pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time;
strata time sex01 race3 agec5; 
run;

proc phreg data=new.mort_temp_ses_all; 
class sex01(ref='1') race3 agec5 time/param=ref;*repeat for each reference group;
model count*dead(0)= sumtemp wintemp sdsumtemp sdwintemp sumtemp*sex01 pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time;
strata time sex01 race3 agec5; 
run;

proc phreg data=new.mort_temp_ses_all; 
class sex01 race5(ref='1') race3 agec5 time/param=ref; *repeat for each reference group;
model count*dead(0)= sumtemp wintemp sdsumtemp sdwintemp sumtemp*race5 pctbelpovgt64 H076001 pop_den pctgreen06 smkrate2 percent_diab percent_ltpia time;
strata time sex01 race3 agec5; 
run;

data new.mort_temp_ses_all;
set new.mort_temp_ses_all;
pop_den01=1;
if pop_den<47.3 then pop_den01=0; *33.3th percentile of poplation density;
run;

proc phreg data=new.mort_temp_ses_all; 
class sex01 race3 agec5 time pop_den01(ref='1') /param=ref;*repeat for each reference group;
model count*dead(0)= sumtemp wintemp sdsumtemp sdwintemp sumtemp*pop_den01 pctbelpovgt64 H076001 pop_den01 pctgreen06 smkrate2 percent_diab percent_ltpia time;
strata time sex01 race3 agec5; 
run;

proc phreg data=mort; 
class sex01 race3 agec5 time age01(ref='0') /param=ref;*repeat for each reference group;
model count*dead(0)=  zipmean_sum zipmean_sum*age01 zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif pctbelpovgt64 H076001 pop_den  pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

proc phreg data=mort; 
class  race3 agec5 time sex01(ref='1') /param=ref;*repeat for each reference group;
model count*dead(0)= zipmean_sum zipmean_sum*sex01 zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif pctbelpovgt64 H076001 pop_den  pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

proc phreg data=mort; 
class  race3 agec5 time race5(ref='1') /param=ref;*repeat for each reference group;
model count*dead(0)= zipmean_sum zipmean_sum*race5 zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif pctbelpovgt64 H076001 pop_den  pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;

proc phreg data=mort; 
class sex01 race3 agec5 time pop_den01(ref='0') /param=ref;*repeat for each reference group;
model count*dead(0)= zipmean_sum zipmean_sum*pop_den01 zipmean_win zipmean_sdsum zipmean_sdwin sumdif windif sdsumdif sdwindif pctbelpovgt64 H076001 pop_den01  pctgreen06 smkrate2 percent_diab percent_ltpia time/rl;
strata time sex01 race3 agec5; 
run;




